Processing Spatial Transcriptomic Data
Larry Chen
2023-01-18
NOTE: This code utilizes Seurat v4.1.1. You will need to install seurat version 4.1.1 using ‘remotes’. We recommend installing the version in separate location from your regular seurat package.
remotes::install_version('Seurat', version = '4.1.1', lib = "C:/Users/Ji Lab/AppData/Local/R/alt_packages/Seurat 4.1.1" )Import Libraries
library(Seurat, lib.loc = "C:/Users/Ji Lab/AppData/Local/R/alt_packages/Seurat 4.1.1")
library(dplyr)
library(ggplot2)sessioninfo::session_info()%>%
details::details(
summary = 'Current session info',
open = TRUE
)Current session info
─ Session info ───────────────────────────────────────────────────────────────
setting value
version R version 4.2.2 (2022-10-31 ucrt)
os Windows 10 x64 (build 22621)
system x86_64, mingw32
ui RTerm
language (EN)
collate English_United States.utf8
ctype English_United States.utf8
tz America/New_York
date 2023-01-19
pandoc 2.19.2 @ C:/Program Files/RStudio/resources/app/bin/quarto/bin/tools/ (via rmarkdown)
─ Packages ───────────────────────────────────────────────────────────────────
! package * version date (UTC) lib source
abind 1.4-5 2016-07-21 [1] CRAN (R 4.2.0)
assertthat 0.2.1 2019-03-21 [1] CRAN (R 4.2.1)
bslib 0.4.2 2022-12-16 [1] CRAN (R 4.2.2)
cachem 1.0.6 2021-08-19 [1] CRAN (R 4.2.1)
cli 3.4.1 2022-09-23 [1] CRAN (R 4.2.1)
clipr 0.8.0 2022-02-22 [1] CRAN (R 4.2.1)
cluster 2.1.4 2022-08-22 [2] CRAN (R 4.2.2)
codetools 0.2-18 2020-11-04 [2] CRAN (R 4.2.2)
colorspace 2.0-3 2022-02-21 [1] CRAN (R 4.2.1)
cowplot 1.1.1 2020-12-30 [1] CRAN (R 4.2.1)
data.table 1.14.6 2022-11-16 [1] CRAN (R 4.2.2)
DBI 1.1.3 2022-06-18 [1] CRAN (R 4.2.1)
deldir 1.0-6 2021-10-23 [1] CRAN (R 4.2.0)
desc 1.4.2 2022-09-08 [1] CRAN (R 4.2.1)
details 0.3.0 2022-03-27 [1] CRAN (R 4.2.2)
digest 0.6.30 2022-10-18 [1] CRAN (R 4.2.1)
dplyr * 1.0.10 2022-09-01 [1] CRAN (R 4.2.1)
ellipsis 0.3.2 2021-04-29 [1] CRAN (R 4.2.1)
evaluate 0.19 2022-12-13 [1] CRAN (R 4.2.2)
fansi 1.0.3 2022-03-24 [1] CRAN (R 4.2.1)
fastmap 1.1.0 2021-01-25 [1] CRAN (R 4.2.1)
fitdistrplus 1.1-8 2022-03-10 [1] CRAN (R 4.2.1)
future 1.30.0 2022-12-16 [1] CRAN (R 4.2.2)
future.apply 1.10.0 2022-11-05 [1] CRAN (R 4.2.2)
generics 0.1.3 2022-07-05 [1] CRAN (R 4.2.1)
ggplot2 * 3.4.0 2022-11-04 [1] CRAN (R 4.2.2)
ggrepel 0.9.2 2022-11-06 [1] CRAN (R 4.2.2)
ggridges 0.5.4 2022-09-26 [1] CRAN (R 4.2.1)
globals 0.16.2 2022-11-21 [1] CRAN (R 4.2.2)
glue 1.6.2 2022-02-24 [1] CRAN (R 4.2.1)
goftest 1.2-3 2021-10-07 [1] CRAN (R 4.2.0)
gridExtra 2.3 2017-09-09 [1] CRAN (R 4.2.1)
gtable 0.3.1 2022-09-01 [1] CRAN (R 4.2.1)
htmltools 0.5.4 2022-12-07 [1] CRAN (R 4.2.2)
htmlwidgets 1.6.0 2022-12-15 [1] CRAN (R 4.2.2)
httpuv 1.6.6 2022-09-08 [1] CRAN (R 4.2.2)
httr 1.4.4 2022-08-17 [1] CRAN (R 4.2.1)
ica 1.0-3 2022-07-08 [1] CRAN (R 4.2.1)
igraph 1.3.5 2022-09-22 [1] CRAN (R 4.2.1)
irlba 2.3.5.1 2022-10-03 [1] CRAN (R 4.2.2)
jquerylib 0.1.4 2021-04-26 [1] CRAN (R 4.2.1)
jsonlite 1.8.3 2022-10-21 [1] CRAN (R 4.2.2)
KernSmooth 2.23-20 2021-05-03 [2] CRAN (R 4.2.2)
knitr 1.41 2022-11-18 [1] CRAN (R 4.2.2)
later 1.3.0 2021-08-18 [1] CRAN (R 4.2.1)
lattice 0.20-45 2021-09-22 [2] CRAN (R 4.2.2)
lazyeval 0.2.2 2019-03-15 [1] CRAN (R 4.2.1)
leiden 0.4.3 2022-09-10 [1] CRAN (R 4.2.1)
lifecycle 1.0.3 2022-10-07 [1] CRAN (R 4.2.1)
listenv 0.9.0 2022-12-16 [1] CRAN (R 4.2.2)
lmtest 0.9-40 2022-03-21 [1] CRAN (R 4.2.1)
magrittr 2.0.3 2022-03-30 [1] CRAN (R 4.2.1)
MASS 7.3-58.1 2022-08-03 [2] CRAN (R 4.2.2)
Matrix 1.5-4 2022-11-14 [1] R-Forge (R 4.2.2)
matrixStats 0.63.0 2022-11-18 [1] CRAN (R 4.2.2)
mgcv 1.8-41 2022-10-21 [2] CRAN (R 4.2.2)
mime 0.12 2021-09-28 [1] CRAN (R 4.2.0)
miniUI 0.1.1.1 2018-05-18 [1] CRAN (R 4.2.1)
munsell 0.5.0 2018-06-12 [1] CRAN (R 4.2.1)
nlme 3.1-160 2022-10-10 [2] CRAN (R 4.2.2)
parallelly 1.33.0 2022-12-14 [1] CRAN (R 4.2.2)
patchwork 1.1.2 2022-08-19 [1] CRAN (R 4.2.1)
pbapply 1.6-0 2022-11-16 [1] CRAN (R 4.2.1)
pillar 1.8.1 2022-08-19 [1] CRAN (R 4.2.1)
pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 4.2.1)
plotly 4.10.1 2022-11-07 [1] CRAN (R 4.2.2)
plyr 1.8.8 2022-11-11 [1] CRAN (R 4.2.2)
png 0.1-8 2022-11-29 [1] CRAN (R 4.2.2)
polyclip 1.10-4 2022-10-20 [1] CRAN (R 4.2.1)
progressr 0.12.0 2022-12-13 [1] CRAN (R 4.2.2)
promises 1.2.0.1 2021-02-11 [1] CRAN (R 4.2.1)
purrr 0.3.5 2022-10-06 [1] CRAN (R 4.2.1)
R6 2.5.1 2021-08-19 [1] CRAN (R 4.2.1)
RANN 2.6.1 2019-01-08 [1] CRAN (R 4.2.1)
RColorBrewer 1.1-3 2022-04-03 [1] CRAN (R 4.2.0)
Rcpp 1.0.9 2022-07-08 [1] CRAN (R 4.2.1)
RcppAnnoy 0.0.20 2022-10-27 [1] CRAN (R 4.2.2)
reshape2 1.4.4 2020-04-09 [1] CRAN (R 4.2.1)
reticulate 1.26-9000 2022-11-29 [1] Github (rstudio/reticulate@a1d7f7f)
rlang 1.0.6 2022-09-24 [1] CRAN (R 4.2.1)
rmarkdown 2.19 2022-12-15 [1] CRAN (R 4.2.2)
ROCR 1.0-11 2020-05-02 [1] CRAN (R 4.2.1)
rpart 4.1.19 2022-10-21 [2] CRAN (R 4.2.2)
rprojroot 2.0.3 2022-04-02 [1] CRAN (R 4.2.1)
rstudioapi 0.14 2022-08-22 [1] CRAN (R 4.2.1)
Rtsne 0.16 2022-04-17 [1] CRAN (R 4.2.1)
sass 0.4.4 2022-11-24 [1] CRAN (R 4.2.2)
scales 1.2.1 2022-08-20 [1] CRAN (R 4.2.1)
scattermore 0.8 2022-02-14 [1] CRAN (R 4.2.1)
sctransform 0.3.5 2022-09-21 [1] CRAN (R 4.2.2)
sessioninfo 1.2.2 2021-12-06 [1] CRAN (R 4.2.1)
VP Seurat * 4.1.1 2022-11-18 [?] CRAN (R 4.2.2) (on disk 4.3.0)
SeuratObject * 4.1.3 2022-11-07 [1] CRAN (R 4.2.2)
shiny 1.7.4 2022-12-15 [1] CRAN (R 4.2.2)
sp 1.5-1 2022-11-07 [1] CRAN (R 4.2.2)
spatstat.core 2.4-4 2022-05-18 [1] CRAN (R 4.2.1)
spatstat.data 3.0-0 2022-10-21 [1] CRAN (R 4.2.2)
spatstat.geom 3.0-3 2022-10-25 [1] CRAN (R 4.2.2)
spatstat.random 3.0-1 2022-11-03 [1] CRAN (R 4.2.2)
spatstat.sparse 3.0-0 2022-10-21 [1] CRAN (R 4.2.2)
spatstat.utils 3.0-1 2022-10-19 [1] CRAN (R 4.2.2)
stringi 1.7.8 2022-07-11 [1] CRAN (R 4.2.1)
stringr 1.5.0 2022-12-02 [1] CRAN (R 4.2.2)
survival 3.4-0 2022-08-09 [2] CRAN (R 4.2.2)
tensor 1.5 2012-05-05 [1] CRAN (R 4.2.0)
tibble 3.1.8 2022-07-22 [1] CRAN (R 4.2.1)
tidyr 1.2.1 2022-09-08 [1] CRAN (R 4.2.1)
tidyselect 1.2.0 2022-10-10 [1] CRAN (R 4.2.1)
utf8 1.2.2 2021-07-24 [1] CRAN (R 4.2.1)
uwot 0.1.14 2022-08-22 [1] CRAN (R 4.2.1)
vctrs 0.5.1 2022-11-16 [1] CRAN (R 4.2.2)
viridisLite 0.4.1 2022-08-22 [1] CRAN (R 4.2.1)
withr 2.5.0 2022-03-03 [1] CRAN (R 4.2.1)
xfun 0.35 2022-11-16 [1] CRAN (R 4.2.2)
xml2 1.3.3 2021-11-30 [1] CRAN (R 4.2.1)
xtable 1.8-4 2019-04-21 [1] CRAN (R 4.2.1)
yaml 2.3.6 2022-10-18 [1] CRAN (R 4.2.1)
zoo 1.8-11 2022-09-17 [1] CRAN (R 4.2.2)
[1] C:/Users/Ji Lab/AppData/Local/R/win-library/4.2
[2] C:/Program Files/R/R-4.2.2/library
V ── Loaded and on-disk version mismatch.
P ── Loaded and on-disk path mismatch.
──────────────────────────────────────────────────────────────────────────────
Load Objects
dir = "C:/Users/Ji Lab/Documents/JID manuscript/andrew_scripts/orig_obj/ST/"
vis.dir = paste(dir,"V10F24-007_A1/outs/",sep="")
vis.dir2 = paste(dir,"V10F24-007_B1/outs/",sep="")
n23v1 = Load10X_Spatial(data.dir = vis.dir, slice = "rep1")## as(<dgTMatrix>, "dgCMatrix") is deprecated since Matrix 1.5-0; do as(., "CsparseMatrix") insteadn23v2 = Load10X_Spatial(data.dir = vis.dir2, slice = "rep2")Merge slices and assign sample identities to cells
n23_p1 = merge(n23v1,n23v2,add.cell.ids = c("rep1","rep2"))
n23_p1$sample = rep(1,ncol(n23_p1))
n23_p1$sample[(ncol(n23v1)+1):(ncol(n23v1)+ncol(n23v2))] = 2cut out bottom portion of spots on folded/noisy tissue
discard_table = read.table(paste(vis.dir2,"Extra.csv",sep=""),sep = ",",row.names = 1,header = T,stringsAsFactors = F)
new_row = paste("rep2_",rownames(discard_table),sep="")
rownames(discard_table) = new_row
discard_spots = new_row[which(discard_table[,"Extra"]=="Discard")]
keep_spots = colnames(n23_p1)[!colnames(n23_p1) %in% discard_spots]
n23_p1 = subset(n23_p1, cells = keep_spots)
n23_p1 = subset(n23_p1, nCount_Spatial > 200)Normalize data
n23_p1 <- SCTransform(n23_p1, assay = "Spatial", return.only.var.genes = F)## Calculating cell attributes from input UMI matrix: log_umi## Variance stabilizing transformation of count matrix of size 17532 by 4141## Model formula is y ~ log_umi## Get Negative Binomial regression parameters per gene## Using 2000 genes, 4141 cells##
|
| | 0%
|
|================== | 25%
|
|=================================== | 50%
|
|==================================================== | 75%
|
|======================================================================| 100%## Found 62 outliers - those will be ignored in fitting/regularization step## Second step: Get residuals using fitted parameters for 17532 genes##
|
| | 0%
|
|== | 3%
|
|==== | 6%
|
|====== | 8%
|
|======== | 11%
|
|========== | 14%
|
|============ | 17%
|
|============== | 19%
|
|================ | 22%
|
|================== | 25%
|
|=================== | 28%
|
|===================== | 31%
|
|======================= | 33%
|
|========================= | 36%
|
|=========================== | 39%
|
|============================= | 42%
|
|=============================== | 44%
|
|================================= | 47%
|
|=================================== | 50%
|
|===================================== | 53%
|
|======================================= | 56%
|
|========================================= | 58%
|
|=========================================== | 61%
|
|============================================= | 64%
|
|=============================================== | 67%
|
|================================================= | 69%
|
|=================================================== | 72%
|
|==================================================== | 75%
|
|====================================================== | 78%
|
|======================================================== | 81%
|
|========================================================== | 83%
|
|============================================================ | 86%
|
|============================================================== | 89%
|
|================================================================ | 92%
|
|================================================================== | 94%
|
|==================================================================== | 97%
|
|======================================================================| 100%## Computing corrected count matrix for 17532 genes##
|
| | 0%
|
|== | 3%
|
|==== | 6%
|
|====== | 8%
|
|======== | 11%
|
|========== | 14%
|
|============ | 17%
|
|============== | 19%
|
|================ | 22%
|
|================== | 25%
|
|=================== | 28%
|
|===================== | 31%
|
|======================= | 33%
|
|========================= | 36%
|
|=========================== | 39%
|
|============================= | 42%
|
|=============================== | 44%
|
|================================= | 47%
|
|=================================== | 50%
|
|===================================== | 53%
|
|======================================= | 56%
|
|========================================= | 58%
|
|=========================================== | 61%
|
|============================================= | 64%
|
|=============================================== | 67%
|
|================================================= | 69%
|
|=================================================== | 72%
|
|==================================================== | 75%
|
|====================================================== | 78%
|
|======================================================== | 81%
|
|========================================================== | 83%
|
|============================================================ | 86%
|
|============================================================== | 89%
|
|================================================================ | 92%
|
|================================================================== | 94%
|
|==================================================================== | 97%
|
|======================================================================| 100%## Calculating gene attributes## Wall clock passed: Time difference of 46.28732 secs## Determine variable features## Place corrected count matrix in counts slot## Centering data matrix## Set default assay to SCTPCA/UMAP/clustering
keep.dims <- 1:20
n23_p1 <- RunPCA(n23_p1, verbose = FALSE)
ProjectDim(n23_p1, reduction = "pca", dims = 1:20)## PC_ 1
## Positive: TCHH, KRT81, S100A3, KRT71, KRT86, KRT85, KRT25, KRT31, PRR9, KRT27
## PSORS1C2, KRT33B, KRT83, KRTAP1-5, KRTAP2-2, KRTAP11-1, KRTAP3-3, KRTAP1-3, KRTAP3-1, LY6G6F-LY6G6D
## Negative: FADS2, KRT79, ALOX15B, MGST1, THRSP, GAL, FASN, AWAT2, HMGCS1, DGAT2
## ACSL1, FAR2, CIDEA, ELOVL5, MSMO1, APOC1, FDPS, FADS1, SOAT1, PM20D1
## PC_ 2
## Positive: DCD, MUCL1, SCGB2A2, PIP, SCGB1D2, SCGB1B2P, LTF, ZG16B, SAA2, SLC12A2
## SAA1, SNORC, KRT7, KRT19, S100A8, AQP5, SFRP1, S100A9, STAC2, S100A7
## Negative: FADS2, KRT79, MGST1, THRSP, ALOX15B, GAL, HMGCS1, DGAT2, AWAT2, FASN
## ACSL1, ELOVL3, FAR2, APOC1, CIDEA, MSMO1, INSIG1, TMEM91, PM20D1, FDPS
## PC_ 3
## Positive: S100A8, S100A9, S100A7, KRT1, KRTDAP, KRT10, CALML5, FABP5, SBSN, CSTA
## FLG, SPRR1B, DMKN, AQP3, LOR, CNFN, SPRR2E, KRT6A, HOPX, LY6D
## Negative: DCD, MUCL1, SCGB2A2, PIP, SAA1, SCGB1B2P, SCGB1D2, SAA2, ZG16B, KRT7
## SLC12A2, SNORC, LTF, KRT19, AZGP1, AQP5, SFRP1, STAC2, CLDN10, KRT8
## PC_ 4
## Positive: DCD, MUCL1, PIP, SCGB2A2, SCGB1B2P, SCGB1D2, ZG16B, SAA1, SAA2, S100A8
## SNORC, KRT7, AZGP1, S100A9, S100A7, LTF, AQP5, KRTDAP, KRT19, MT-CO1
## Negative: CCDC80, C3, IGFBP5, IGKC, IGFBP4, CXCL14, DCN, IGHG1, COL1A1, COL1A2
## C1R, MAP1B, COL6A2, APOD, EBF1, SPARCL1, COL6A3, IGHG3, CXCL12, NEAT1
## PC_ 5
## Positive: S100A2, KRT14, KRT5, KRT6B, TMSB4X, KRT17, PDZRN3, KRT75, KRT6A, GJA1
## KRT16, LGALS1, S100A14, S100P, S100A10, SERPINA1, CTSV, TMSB10, MRC2, MT2A
## Negative: FLG, FLG2, LOR, CALML5, KRT1, CNFN, HOPX, SPRR2E, SLURP1, KRT10
## ARG1, KRTDAP, CSTA, ASPRV1, DCD, SBSN, LCE3D, HAL, IGKC, KRT2
## PC_ 6
## Positive: KRT85, TCHH, KRT35, MT4, KRT71, KRT32, KRT25, KRT31, FABP4, KRT82
## KRT27, SLPI, KRT73, S100A3, LY6G6F-LY6G6D, PRR9, SP6, LAP3, KRT28, SELENBP1
## Negative: KRT6B, KRT14, KRTAP4-12, KRTAP4-7, KRT5, KRTAP2-2, KRTAP4-9, KRTAP4-6, KRTAP3-2, KRTAP9-4
## KRTAP1-3, KRTAP1-1, KRTAP3-1, KRTAP9-3, KRTAP4-2, KRTAP2-1, KRTAP4-3, KRTAP4-11, KRTAP4-4, PDZRN3
## PC_ 7
## Positive: PDZRN3, CTSV, KRT6B, MUCL1, KRT75, IL18, TPM1, KLK5, KLK7, SERPINA1
## CST6, SPINK5, S100P, ASPRV1, MYL9, CRCT1, LGALS1, GJA1, CDSN, LY6G6C
## Negative: AQP3, S100A8, LY6D, IGKC, IGHG1, S100A9, TYRP1, FGFBP1, MT2A, IGHG3
## COL17A1, MT1E, DCT, PMEL, IGHG2, IGHG4, KRT10, IGLC3, IGHA1, B2M
## PC_ 8
## Positive: MT2A, IGFBP5, DCN, COL1A2, COL1A1, CTSL, S100A8, C1R, CCDC80, FOS
## AQP3, CFH, S100A9, IFITM3, SPARC, COL3A1, FSTL1, COL6A2, TIMP1, LY6D
## Negative: IGHG1, IGKC, IGHG3, IGHG2, IGLC3, IGHA1, IGHG4, IGLC2, MZB1, JCHAIN
## DCD, CST6, SPINK5, CD79A, CTSV, KRT6B, DERL3, IGHGP, FLG, MUCL1
## PC_ 9
## Positive: TCHH, FLG, SPINK5, NEAT1, DSP, DCD, HSP90AA1, CALD1, RBM25, KRT85
## KRT2, PRRC2C, KRT10, KRT1, EIF5B, KRT81, MYH11, LUC7L3, PABPC1, DST
## Negative: IGKC, IGHG1, IGHG3, IGHG2, C3, TMSB4X, IGFBP4, CD74, CXCL12, PSAPL1
## IGLC2, HLA-DRA, DCN, VSIG8, JCHAIN, C1S, APOD, SPRR2D, IGLC3, VIM
## PC_ 10
## Positive: COMP, CSPG4, BASP1, SLC7A8, LGALS1, COL16A1, TNC, KRT74, LRRC15, THBS1
## CDSN, BGN, SELENOP, KRTAP10-10, LGR5, APMAP, KRTAP10-1, CREB5, TIMP3, FGF18
## Negative: CST6, PSAPL1, KRT17, KRT85, KRT6B, MUCL1, KRT81, SEC14L6, KRT31, CALML3
## KRT6A, GAL, ELOVL4, KRT16, SPINK5, CTSV, MUC1, CRAT, KRT86, SOAT1
## PC_ 11
## Positive: TCHH, KRT25, KRT71, KRT35, KRT27, FABP9, KRTAP9-3, C3, SLPI, KRT85
## TCHHL1, KRTAP4-6, KRTAP4-9, CTSL, CXCL14, KRTAP4-3, KRTAP4-8, KRTAP1-5, KRTAP9-4, KRTAP4-12
## Negative: DES, MYH11, MYL9, CNN1, TPM2, TAGLN, SYNPO2, MAP1B, ACTG2, KRTAP11-1
## CAPN12, OXCT2, SORBS1, ACTA2, ENGASE, FLNC, CST6, ACTA1, CALD1, CALHM4
## PC_ 12
## Positive: CXCL14, FLG, COL1A1, COL1A2, THBS2, SEC14L6, KRT1, KRT15, LTF, PSAPL1
## KRT10, COMP, KRTAP17-1, CREB5, ELOVL4, KRT5, KRT81, KRT83, AADACL3, SFRP1
## Negative: MYL9, DES, MYH11, TAGLN, KRT35, CNN1, TPM2, ACTA2, C11orf96, ACTG2
## SYNPO2, FLNC, ACTA1, MYLK, MT4, S100P, KRT85, IGHG1, ITGA5, SPARCL1
## PC_ 13
## Positive: TCHH, PSAPL1, GAL, SEC14L6, ELOVL4, MUC1, CRAT, AADACL3, ELOVL1, KRT71
## SOAT1, S100A9, COMP, FABP9, KRT25, MYL9, KRT27, IGHG4, S100A8, DES
## Negative: APMAP, MUCL1, S100P, AWAT2, ACO1, TF, PDZRN3, FDPS, CTSV, HSD3B1
## FASN, CYP4F8, CST1, KRT16, SRD5A1, GLDC, PGRMC1, C3, EFHD1, APOD
## PC_ 14
## Positive: KRT35, MT4, PMEL, KRT15, PSAPL1, SEC14L6, APCDD1, GPNMB, ELOVL4, GAL
## AADACL3, SOAT1, ID3, SELENBP1, KRTAP4-11, LEF1, FLG2, KRTAP4-7, DAPL1, MLANA
## Negative: TCHH, KRT81, KRT17, S100A9, S100A8, S100A7, KRT6A, SPRR2A, KRT6C, S100A7A
## KRT83, PI3, KRT6B, SERPINB4, KRT33B, KRT16, KRTAP11-1, S100A2, APMAP, AWAT2
## PC_ 15
## Positive: KRT10, KRT1, FABP5, KRT2, KRTDAP, LGALS7B, KRT16, LY6D, SLPI, PDZRN3
## SERPINA1, KRT6C, C3, KRT75, CSTA, IGFL1, SPRR1A, TPM1, S100A8, MYL9
## Negative: CST6, MT2A, SPRR4, MMP3, KRT17, SPRR2E, SPRR2G, LCE3D, ATP6V1C2, IER3
## SPINK5, SPRR2B, KLK6, PTN, KRT77, S100A2, LCE3E, PRSS22, IGFL2, C15orf48
## PC_ 16
## Positive: DCD, KRT74, KRTAP10-10, NEAT1, CAPN8, DUSP5, KRT6A, LRRC15, KRTAP10-1, MUCL1
## LGALS3, KRT72, KRT27, KRT25, PINLYP, KRTAP10-3, LYPD5, SAMHD1, KRT26, SPRR2D
## Negative: KRT81, KRT31, KRT86, KRT83, COMP, KRTAP11-1, KRT33B, TNC, LY6G6F-LY6G6D, KRTAP2-1
## CSPG4, THBS1, THBS2, KRT85, CREB5, RNF152, S100A2, SLURP1, KRTAP1-3, KRTAP8-1
## PC_ 17
## Positive: CLDN5, PECAM1, A2M, VWF, AQP1, TCHH, SOX18, CD93, ACKR1, SPARCL1
## EGFL7, RGS5, LTF, NR2F2, ADGRF5, CD74, TGFBR2, RNASE1, ESAM, EPAS1
## Negative: DCD, COL1A2, COL1A1, CCDC80, IGKC, DCN, CXCL14, IGLC3, FBN1, IGHG1
## SCGB2A2, FBLN1, MUCL1, IGHG4, PIP, TNXB, AEBP1, DES, IGHA1, THBS1
## PC_ 18
## Positive: DCD, MUCL1, SCGB2A2, PI3, KRT6C, PIP, S100A2, SCGB1D2, S100P, SPRR2A
## TNC, TIMP3, THBS1, KRT85, ANGPTL7, FGF18, GJB6, FOS, MMP7, COMP
## Negative: CST6, SPRR4, TCHH, SPINK5, LCN2, C15orf48, CTSV, KRT15, IGFL2, FABP9
## CLU, KRT17, LTF, SAA1, CALHM5, SAA2, RPL28, KRT19, RPS18, PTN
## PC_ 19
## Positive: LTF, S100A2, TYRP1, MMP3, COL17A1, PI3, AQP3, MT1E, DCT, ITGA6
## S100A6, MT2A, KRT14, SAA1, SPRR2G, KRT6C, IGKC, PHLDA1, KRT16, SAA2
## Negative: CST6, KRTDAP, SPRR4, DCD, COMP, FABP5, S100A7, SCGB1D2, SCGB2A2, KRT77
## THBS1, PIP, KRT10, CLDN5, TSC22D1, CALML5, IGFL2, DEFB1, LGALS7B, CSPG4
## PC_ 20
## Positive: KRT35, MT4, PI3, S100A7A, S100A8, S100A7, KRT17, SPRR2A, KRT85, SERPINB4
## SPRR1B, SELENBP1, COMP, KRT6C, SAA1, STAC2, S100A9, DEFB4A, NDRG1, SAA2
## Negative: TCHH, COL17A1, TYRP1, DCT, MUCL1, IL18, FLG2, IFITM3, KRT14, SPINK5
## MT1E, KRT5, FGFBP1, FLG, KRT81, KRT15, AQP3, KRT33B, DCD, FABP9## An object of class Seurat
## 51070 features across 4141 samples within 2 assays
## Active assay: SCT (17532 features, 3000 variable features)
## 1 other assay present: Spatial
## 1 dimensional reduction calculated: pca
## 2 images present: rep1, rep2n23_p1 <- FindNeighbors(object = n23_p1, dims = keep.dims, verbose = FALSE, reduction = "pca") #PCA
n23_p1 <- FindClusters(object = n23_p1, verbose = FALSE)
n23_p1 <- FindClusters(object = n23_p1, resolution = 1, verbose = FALSE)
n23_p1 <- RunUMAP(object = n23_p1, dims = keep.dims, verbose = FALSE, reduction = "pca")## Warning: The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
## To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
## This message will be shown once per sessionGraphing
DimPlot(n23_p1,label=T)
SpatialDimPlot(n23_p1, group.by = "SCT_snn_res.0.8", label = TRUE, label.size = 3)## Scale for fill is already present.
## Adding another scale for fill, which will replace the existing scale.
## Scale for fill is already present.
## Adding another scale for fill, which will replace the existing scale.
Idents(n23_p1) = 'SCT_snn_res.0.8'
SpatialDimPlot(n23_p1, label = TRUE, label.size = 3)## Scale for fill is already present.
## Adding another scale for fill, which will replace the existing scale.
## Scale for fill is already present.
## Adding another scale for fill, which will replace the existing scale.
DimPlot(n23_p1,group.by = "sample")
VlnPlot(n23_p1,features = c("KRT10","KRT1"))
Note: clustering may differ slightly from what’s shown in figure 2d, due
to instability in clustering / differences in seurat/ seurat dependency
versions
SpatialDimPlot(n23_p1, images = "rep2")
FeaturePlot(n23_p1,features = c("COL17A1","SPRR2D","PTHLH","MMP10"))
VlnPlot(n23_p1,features = c("MMP10","COL17A1","IL24","KRT15","PTHLH","CCL2"))
res = "SCT_snn_res.0.8"
res = "SCT_snn_res.1"
Idents(n23_p1) = res
n23_p1.markers <- FindAllMarkers(object = n23_p1, only.pos = T)
n23_p1_pca_noreg.markers = n23_p1.markers
n23_p1.markers_sig = subset(n23_p1.markers, p_val_adj<0.05)
n23_p1.markers_res1 <- FindAllMarkers(object = n23_p1, only.pos = T)
n23_p1.markers_res1_sig = subset(n23_p1.markers_res1, p_val_adj<0.05)setwd("~/Dropbox/NS_scRNA-seq/ST/N23_seurat/n23_p1")
write.table(n23_p1.markers,file="n23_p1_res0.8_seurat_markers.csv",sep=",",row.names = T, col.names = T)
write.table(n23_p1.markers_res1,file="n23_p1_res1_seurat_markers.csv",sep=",",row.names = T, col.names = T)
write.table(n23_p1.markers,file="n23_p1_pca_no_reg_res0.8_seurat_markers.csv",sep=",",row.names = T, col.names = T)n23_p1.pca.markers <- FindAllMarkers(object = n23_p1, only.pos = T)
write.table(n23_p1.pca.markers,file="n23_p1_pca_res0.8_seurat_markers.csv",sep=",",row.names = T, col.names = T)
n23_p1.markers = read.table("n23_p1_res0.8_seurat_markers.csv", sep = ",",header=T,row.names = 1, stringsAsFactors = F)
n23_p1.markers_sig = subset(n23_p1.markers, p_val_adj<0.05)
match(tme_tsk_ligs,n23_p1.markers_sig$gene)
tme_tsk_ligs[tme_tsk_ligs %in% n23_p1.markers_sig$gene]
n23_p1_top10 = n23_p1.markers %>% group_by(cluster) %>% top_n(10, avg_logFC)
DoHeatmap(n23_p1,features = n23_p1_top10$gene)
n23_p1_pca_top10 = n23_p1.pca.markers %>% group_by(cluster) %>% top_n(10, avg_logFC)#PCA
setwd("~/Dropbox/NS_scRNA-seq/ST/N23_seurat/n23_p1/pca/")
clus = levels(Idents(n23_p1))
for (i in 1:length(clus)) {
p1 = SpatialDimPlot(n23_p1,cells.highlight = WhichCells(n23_p1,idents=clus[i]))
png(paste("n23_p1_SCT_pca_res0.8_highlight_clus",clus[i],".png",sep=""),width = 8, height=4, units = "in", res = 300)
print(p1)
dev.off()
}Feature plots
SpatialFeaturePlot(n23_p1, features = c('KRT5'), images = 'rep2')
SpatialFeaturePlot(n23_p1, features = c('KRT17'), images = 'rep2')
SpatialFeaturePlot(n23_p1, features = c('MGST1'), images = 'rep2')
SpatialFeaturePlot(n23_p1, features = c('DCD'), images = 'rep2')
SpatialFeaturePlot(n23_p1, features = c('KRT10'), images = 'rep2')
SpatialFeaturePlot(n23_p1, features = c('KRT25'), images = 'rep2')
SpatialFeaturePlot(n23_p1, features = c('AWAT2'), images = 'rep2')
SpatialFeaturePlot(n23_p1, features = c('PIP'), images = 'rep2')
SpatialFeaturePlot(n23_p1, features = c('KRT5'), images = 'rep2')
save object
setwd("C:/Users/Ji Lab/Documents/JID manuscript/andrew_scripts/orig_obj/")
saveRDS(n23_p1,file="n23_p1.Rds")